This kind of stuff really makes me appreciate black's approach to lists as "do what'll minimise git diffs" haha.

This is pre-existing, but global_phase is not the global phase until its abs happens to be 1. This test is actually testing whether the gate is close to a scaled identity. Same with diag_is_identity_up_to_global_phase below.

Maybe global_phase should be renamed global_scale if it can't be constrained always to have magnitude 1. And likewise, this function should be renamed. But maybe these can wait for a rewrite.

A bit inconvenient that the product of "empty iterator" doesn't produce the unit tuple in this form - I think Python's itertools feels more mathematically natural here.

If we really wanted to avoid the large code duplication we could pull out the block of the if into a

fn apply<T: IntoIterator<Item = [usize]>>(...) {}

(or whatever the right iterator item is) and call it twice, but I don't really think it's worth bothering.

These numbers are actually just $-i$ and $i$, which makes it all the funnier to me that the Python-space version defined a full "get the RZ matrix" method lol

Lol, oh I didn't even see that I just was in mechanical porting mode. Do you think it's worth making that explicit here, or just leave it for LLVM to optimize?

I'm fine leaving all the little things like this - larger-scale real performance improvements in this would come from a more complete refactor of the implementation as a whole, and as you say, the compiler can probably constant-fold this anyway.

I'd prefer that this were somehow called the 2-norm or p-norm for p=2 ("somehow" with a valid identifier). l2 is a space, usually a sequence space. But upon googling for uses of p norm, l2 norm, etc. it seems that, as far as the internet is concerned, especially machine-learning internet, lp, Lp, can mean whatever you want them to mean. [EDIT. looks like this is translated from np.linalg.norm, which, along with Julia's LinearAlgebra, calls this the 2-norm)]

We might want a function for p norms with p as a parameter, with 2 as the default, and make sure constants are propagated, or whatever it takes to get the correct norm at compile time. But, unless we have an existing use for p != 2 that can be postponed indefinitely.

In any case, at some point, we will wish that functions like this had been collected in a repo-wide module. (probably not in the ten thousandth "utils.xxx"). And in that case it's worth having somewhat more general trait bounds. Not repeating this function could be important for correctness and performance. For example, the Julia 2-norm function appears to sometimes scale the elements to avoid overflow. However, the function mynorm(v) = sqrt(sum(x -> x^2, v)); is more than 3 times faster than top-level entry point norm. We'd want an api and implementation that reflects our common uses, which might be almost all 2- and 4-element vectors. Indeed, an optimization, which we don't need to make at the moment, would take advantage the length of the vector, which is known here at compile time to be two.

The main point for now is to consider at least collecting these somewhere where they are visible so that the questions of performance and correctness can be more easily addressed at some point

I don't want to get into a bikeshedding conversation about naming, but l2_norm or lp_norm is a pretty standard terminology for this. If you look at the comment in the Julia code below what you linked they use it there too:

https://github.com/JuliaLang/julia/blob/68da780b93518204d874410307791702d5200e29/stdlib/LinearAlgebra/src/generic.jl#L496

# Compute L_p norm ‖x‖ₚ = sum(abs(x).^p)^(1/p)

That being said the response to the main comment here is that we can look at creating a dedicated mathematical functions module, or even a separate crate if it's generic enough for these kind of things in a follow up PR. Ideally we'd contribute this to something like ndarray or faer imo though assuming a generic implementation. I only did it like this because it was so small and simple and the usage was very minimal. So far the only function I think that falls into this category is the 2x2 determinant function which is used in this PR, the two qubit decomposer, and the one qubit decomposer (where it was originally added). But again it's just a single line so I didn't feel like it was worth creating a dedicate module for it.